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We find a way to write the perturbation equation in a Kerr background as a coupled system of one 
dimensional equations for the different modes in the time domain. Numerical simulations show that 
the dominant mode in the gravitational response is the one corresponding to the mode of the initial 
perturbation, allowing us to conjecture that the coupling among the modes has a weak influence in 
our one dimensional system of equations. We conclude that by neglecting the coupling terms it can 
be obtained a one dimensional harmonic equation which indeed describes with good accuracy the 
f*"^ ■ gravitational response from the Kerr black hole. This result may help to understand the structure 

CNj ' of test fields in a Kerr background and even to generate accurate waveforms for various cases in an 

< i , efficient manner. 

a 



in 



5-H . 



> 

(N 



X 



PACS numbers: 11.15.Bt, 04.30.-w, 04.25.dg, 02.30.Jr, 95.30.Sf 



I. INTRODUCTION 



There is an intrinsic interest on understanding the structure of black holes and the behavior of the fields close to 
them, where the gravity effects are strong. One natural way to describe such a system is through the linearization of 
the Einstein equations around a fixed background solution. From the theoretical point of view, linear perturbations 
allow us to study the stability of the spacetime, the damping of the perturbations at late times (ie, the Quasi- Normal- 
cy Modes and the tails) and other phenomena like the super radiance or the extraction of rotational energy from the 
black hole. On a more practical side, it can be used as a tool for computing the gravitational wave emission when 
the spacetime is perturbed by a less massive source. In these cases the spacetime is dominated by the black hole and 
£vq [ it can be considered as a fixed background. A very likely occurrence example is the capture of compact stellar-sized 
£S) ■ objects (ie, a black hole or a star) by the super massive black holes lying at the center of a galaxy. This type of 
Extreme Mass Ratio Inspiral (EMRI) events will generate gravitational waves within the sensitivity band of LISA, 
. the prospected space-based GW detector. In spite of the recent breakthroughs done by Numerical Relativity in the 
f*^ ' full 3D evolution of binary black holes, this is one of the problems which seems not possible to be solved easily by 
t-H . these methods, due mainly to the very different scales presents on the problem. 

One of the frames to study the perturbations is by means the linearization of the Einstein Equations through the 
perturbed Weyl scalar vf^ The perturbations of the Schwarzschild solution (i.e., non-rotating BH) are well known 
since the seventies 0]. In this case, the solution of the angular part is simply given by the spin weighted spherical 
harmonics, while the time-radial part describes the evolution of multipoles which evolve independently. 

Perturbations around Kerr spacetime (i.e., rotating BH) arc much more complicated, since the spacetime is not 
spherically symmetric and the fields can not be decomposed into independently evolving multipoles in the time 
domain. Already Teukolsky [l[ , with a null tetrad of the Kinnersley type, arrived to a separable master perturbation 
equation for the Kerr black hole in Boycr-Lindquist coordinates. However, the angular part was described by the so 
called spheroidal spin weighted harmonics, for which there is no known analytical expression. Moreover, the multipolc 
decomposition could only be performed in the frequency domain, but with separation constants which were functions 
of frequency, implying a non-trivial (and unknown) multipole interaction in the time domain. A expansion on the 
spin parameter was used in Q to separate the modes and study the effects of the mode coupling on the tails, the very 
late behaviour of the perturbations after the ring down. However, their results are by construction restricted to small 
deviations from Schwarzschild spacetime (see also [![)■ 

In this work, we describe the Kerr black hole in horizon penetrating coordinates and define a new tetrad such that 
the main angular operator, in the evolution equation for the Newman-Penrose scalar ^4, has as eigenfunctions the 
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usual spin weighted spherical harmonics. This allows us to write the ^4 in a base of these harmonics, without the 
need to define the spheroidal ones. The action of the main angular operator on ^4 can be expressed in terms of 
well known eigenvalues. This procedure leads to an operator with only radial and temporal derivatives, but still with 
angular coefficients. By using the normalization properties of the spherical harmonics, it is possible to arrive to a 
purely radial-temporal operator for a combination of modes, reducing in this way the perturbation problem in Kerr 
to a coupled system of one dimensional equations for the different modes. A similar procedure can be performed with 
the Maxwell and the Klein- Gordon equations. 

Furthermore, our procedure allows to get a clearer description of the modes composing the gravitational wave in 
terms of the initial data. We anticipate that the dominant mode in the waveform is the one corresponding to the 
mode of the initial data, implying that the coupling terms are weak with respect to the main harmonic operator. This 
result further reduces the problem of generating wave forms in a perturbed Kerr spacetime, since up to some level of 
accuracy is enough to simply solve the equation for the corresponding mode. 

The work is organized as follows: in section UH we present a derivation of the perturbation equation for the 
perturbed ^4 Weyl scalar for a Kerr black hole described in Kerr-Schild penetrating coordinates, and explicitly 
derived the evolution equation for $ = rVE^. In section ITU the function $ is decomposed in terms of spin weighted 
spherical harmonics. After some manipulations, the perturbation equation can be separated in radial-temporal and 
angular parts, resulting in a coupled system of radial-temporal equations. In section IIVI it is shown a first order 
reduction of the coupled system of equations. It is also explained the standard approach, which consists on using a 
different base of angular functions and solve the Teukolsky equation in 2D. In section [V] we describe the numerical 
code used to evolve the two previous systems and present few numerical examples of the gravitational response of the 
black hole due to a gravitational pcrtubation. Finally in section [VTl we present a discussion of the results obtained. 



II. PERTURBATION EQUATION IN KERR SPACETIME 

The derivation of an evolution equation for the perturbed Weyl scalar requires first the choice of a convenient null 
tetrad Z a ^ . Following a previous work [Bjj, we define a null tetrad with two real future-directed null vectors Z^ L = 
and = fc M , where is pointing outward and fc M pointing inward in an hypcrsurface of constant time. The other 
two null vectors lay in the plane perpendicular to the light cone and will be denoted by a complex vector m M , such 
that Z 2 ^ = m* 1 , Z-^ = m*^. These null vectors must satisfy the following equations: 

Z a M Z bfl = rjab, g^y = 2ri ab Z a ^ Z bu ^, (1) 

being rjai, a matrix where all the coefficients are zero except 7712 = 7721 = — %4 = — V13 = V- The constant 77 is equal 
to one if the spacetime has signature {+, — , — , — }, and equal to minus one for opposite signature { — , +, +, +}. 

Following the procedure and definitions described in |5|, Q, we obtain the following evolution equation for the 
perturbed Weyl scalar ^4^ , for type D spaces in vacuum backgrounds, 

((A + 77 (4/i + //+37-7*)) (D - 77 (p- 4e)) - (<f* + (3a + (3* + 4tt - r*)) (<5 + 77 (4/3 -r)) 
- 3 7 n2 1< 2 ) * 4 (1) =v^tvT (2) 

where we have used geometric units. This evolution equation for the perturbed Weyl scalar was originally derived by 
Teukolsky [l|, modulus the coefficient 77 which takes into account the signature convention. It was successfully used 
to describe, in the frequency domain, the gravitational waves generated by a perturbed rotating black hole. 

Now we proceed to use the perturbation equation ([2]) in a background given by the Kerr black hole described in 
penetrating coordinates. In these coordinates, the line element is given by: 

/ M r\ Mr Mr ( Mr\ ( Mr\ 

ds 2 = - f 1- 2-^- J dt 2 + 4-^-didr- 4 a sin 2 9 dt dip + ( I + 2— j- j dr 2 - 2 ( I + 2-^ j a sin 2 9 dr dip + 

Y,d9 2 + [r 2 + a 2 + 2— a 2 sin 2 9 J sin 2 9 dip 2 (3) 

with M is the mass of the black hole, a is its angular momentum per unit mass, and we have defined E = r 2 + a 2 cos 2 9. 
Notice that we have already chosen our signature so that 77 = — 1 . 

The choice of the null tetrad is ambigous, since it is only restricted by the eqs. Q with 77 = —I. Based on the 
procedures used in 0, 0] for the Schwarzschild case, we define the following tetrad: 

= ^-(r 2 + a 2 + 2 Mr, A,0, 2 a) , = k Q (I, -I, 0, 0) , 777^ = ^ (i a sin (9, 0, 1, i esc 9) , (4) 

2 E y 2 (r — ia co&9) 
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where we have used the shortcut A = r 2 + a 2 — 2 Mr. In the spatial asymptotic region (i.e., r — > oo), the real vectors 
of the null tetrad ((31) take the form expected in a flat spacetime / M = i (1, 1,0,0) , fc M = fc (1,-1,0,0). Regarding 
the complex null vector, our choice differs slightly from the one commonly used in the literature [l|, changing the 
denominator by its conjugate complex. This rotation will have a big impact on simplifying the final perturbation 
equation. 

From the definitions of the spinor coefficients Q it is straighforward to obtained that the un-perturbed Weyl scalars 
are 4*0 = = ^3 = ^4 = and \&2 = M/1 3 . The only non- vanishing spinor coefficients in this case are: 

P = r-ia cos8 > P = TS W € = P 2ST' 

- = ^W ; r = i^^; «=f7§^; (3 = -a*+r. (5) 

A straightforward substitution of these quantities in the perturbation equation leads to the following equation for 
the perturbed scalar of ^4: 

^ + * A W=4vT. (6) 



where 



□* = -(r 2 + a 2 cos 2 + 2Mr) * + A ^ + 4Mr + 2 a ■ 



dt 2 dr 2 dtdr drdip 

2 (2r + 37\/ + 2iacos6») — + 6 (r - M) — +4 (7) 

dt or 

d 2 Id 2 „ d „ cos6 d l + cos 2 6> 

oO sm 2 6 dtp 2 d9 sm 2 9 dip sm 2 9 

Notice how remarkable simple turns out to be the operator acting on Vt^ 1 - 1 , with an angular operator \3$ v which is 
exactly the same as the one for the Schwarzschild case. In what follows, we will consider only vacuum spacetimes 
(i.e., T = 0). The explicit form of the operators acting on the source term arepresented in the appendix fA"l 

Our operator (J6]{8]) has exactly the same form as the operator obtained in [8j. However, in their case the operator 
is acting on the function (r — ia cos 9) 4 due to a different choice of the null tetrad. In our case the operator 

is acting directly on ^4^ since our tetrad takes the expected form in the spatial asymptotic region. This implies 
that the Peeling theorem can be applied, recovering the expected decay \&4^ ~ 1/r at far distances. It is then more 
convenient to work with the function $ = r v I / 4^, which will have a constant behavior in the regions far from the 
black hole. The perturbation equation ^ for this quantity takes the form 

-L pt + n e(p ] $ = 0. (9) 

There is only a change in the radial-temporal operator, which now is written as 

d 2 . d 2 .,d 2 . d 2 



nf r = -(r 2 +a 2 cos 2 9 + 2Mr) — + A — + 4Mr — - + 2 a 



dt 2 dr 2 dtdr drdip 

2(2r + M + 2iacos0)|- + 2 ( 2r -M-£) A_ 2 «^ +2 ^±^!. (10 ) 
dt \ r J dr r dip r z 

The main advantage of this perturbation equation with respect to previous ones found in the literature is that the 
angular operator ([5) can be expressed in terms of the cth operators @, US], just as in the Schwarzschild case [f|: 

□^ = 8-1 3-2, (11) 

This implies that the eth operators act initially on a function of spin weight —2. Since the 8 operator raises the spin 
weight and the § one lowers it, their combined action finally give again a function of spin weight —2, so that Y_ 2 l ' m 
is an eigenfunction of such operator: 

Ue v Y-2 hm = B-i 3-2 Y-2' m = -(1-1) (1 + 2) Y. 2 l ' m . (12) 

In addition, we will take advantage from the fact that they arc also eigenfunctions of the azimuthal operator, namely 

^-Y^ m = imY_ 2 l ' m . (13) 
dip 
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These results provide us with a natural way to choose the spin-weighted spherical harmonics as a basis to generate 
the space of solutions for the perturbation equation. This approach is different from the standard one, which started 
by Teukolsky [l[ and have been used by many others authors. In the standard approach, the spin- weighted spheroidal 
harmonic functions are introduced as an eigenfunctions of the angular operator containing the rotation parameter of 
the black hole. The spin- weighted spherical harmonics would correspond to a subset of these functions. The use of 
this basis allows us to know the exact form of the constant that replaces the angular part, so that it is no longer an 
unknown to be determined by the expansions procedures, as for instance in (Tlj . 



III. THE ONE DIMENSIONAL SYSTEM OF TEUKOLSKY EQUATIONS 



In order to separate the angular dependence we first express the perturbation function $ in terms of the standard 
spin- weighted spherical harmonics: 



(14) 



lm 



with / > 2 because the gravitational waves are quadrupole and higher. By using that the spherical harmonics are 
eigenfunction of the angular (|12[) and the azimuthal operator (|13|) we can get almost a radial-temporal equation, 
namely 



2 2 r - M 



Q2 q2 q2 

(r 2 + a 2 cos 2 6 + 2Mr) j— + A—+AMr — — 
v ' at z or z otor 







d 

2 (2r + M + 2ia cos0) — - 



2im- + 2 = (1-1) (1 + 2) 

or r r A 



0. 



(15) 



There are only left two terms with angular dependence, which are proportional to cos 9 and cos 2 9. Following the 
usual procedure to obtain an equation for the coefficients -Rz jln , we multiply the previous equation by the complex 

conjugate Y^ 2 m = Y 2 ~ m , and integrate over the angles: 



J2 j dQY/'- m 'Y_ 



) S+ A ^+ 4M ^ + 2 (2r + M)| 



2 \2r-M 



d_ 

dr 



2im- + 2 

r 



dt 2 dr 2 
AI r + a 2 



(Z-l) (1 + 2) 



dt 



Ri, 



dQY 2 l ''- m 'Y_ 2 l ' m CO s 2 9^R w 

ot z 



dQY 2 1 ' - m ' Y_ 2 l ' m coaO^Rt, 



0. 



(16) 



where we have splitted the sum to isolate the coefficients which include angular functions. 

By using the orthonormalization relation of the spin-weighted spherical harmonics, the angular integral of 



Y 



I 1 , — m' 



Y- 



l,m 



is reduced to Si v 5„ 



On the other hand, the cosinus appearing on the last two terms can be 



written in terms of the spin-weighted spherical harmonics, namely 



cos 2 9 = 



Y 



2,0 



cos9 = 2J-Y 1 ' . 



(17) 



Substituting these relations in the eq. ([To]) lead to two angular integrals with three spin- weighted spherical harmonics. 
These angular integrals can be solved by means of the Wigner 3 — lm symbols (see [l2| and references therein) : 



Y 



h - m Y^ m 'Y*' dn = i, l 



-2 J 
3 fs 



Y'' v dn = j -?-5, 



4^ ^m,m' 



B^m 5l,l'-l + l(l+i) $1,1 



A 



*,./_, + 16- 



l+2,m Ol,l'-2 



(1+2)1 

16 + E~l,m 5l,l' + l + M,m 5l,l'-2 



Bl+l, m 5l,l' + l 

R r | 2(i+4)(Z-3)(Z(Z+l)-3m 2 ) <■ 

^>i+l,m i,I'-l+ 3(2l+3)(l+l)(l)(2l-l) U'~T 



(18) 



In this way, it is obtained a ID coupled system of evolution equations for the gravitational perturbations in Kerr 
spacctimes, which can finally be written as 



-a 2 Ai + 2,, n du Ri+2,m — 4 a Bi + i m ^4 a jj^yidtt — i dtj Ri+i, m + D;,m Ri, 
4aBj, m ^4 a d tt — i OA Ri-i, m — a 2 A^ m d tt Ri~2,m = 0, 



(19) 
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where we have dropped the primes in Z, m and we have defined 



1 / (I - 2)(j - 3)(Z + 2)(Z + l)(j + m)(Z + m - 1)(Z - m)(l - m - 1) 

z,m " I (2/ — 1) (| — 1) V (2Z + 1)(2Z - 3) ' 1 ' 




< (l-2)(l + 2)(l + m)(l-m) 

(2Z + 1)(2Z - 1) ' 1 } 

a 2 ( oG + 4)(Z-3)(Z(Z + l)-3m 2 )\\ 3 2 A d 2 , „ 9 2 

m = "l r +2Mr+ 3 ( 1 + 2 (2I + 3)(I+i)I(2l-l) )) ^ + A ^ +AMr dWr + < 22 > 

, r - m \ 9 „ / a 2 \ a „ a M r + a 2 ,, , ,, . 

M + 4ia 77j — - — + 2 2r-M +iam) — -2im- + 2 5 (Z - 1) (Z + 2) . 

Z(Z + 1) J dt \ r J dr r r 2 v ' y ' 

The angular part dependence has been replaced by a coupling among modes, achieving the goal of reducing, in the 
time domain, the perturbation equation in a Kerr background to a radial-temporal system of (coupled) equations. 

The angular part of the gravitational perturbation is described by the spin-weighted spherical harmonics, while 
the radial-temporal part is obtained by solving the system of equations (|19[) . This system shows the couplings and 
interdependence of the modes. The Ri tTn mode is coupled via the Z number with their neighbors, two above and two 
below. The modes are not coupled with respect to the m— mode. The coupling coefficients A^ m and i?; jTO grow with 
Z, indicating an asymmetry on the coupling. They tend to asymptotic finite values which may be used to bound 
the coefficients, lim Ai :Tn = i and lim Bi tTn — \. They depend symmetrically on m and grow with its modulus. 

Finally, we see that these couplings are determined by the rotational parameter a of the black hole, and vanish in the 
Schwarzschild case. 

Another interesting property of this system is that they are coupled through the temporal derivatives of the modes, 
while the radial dependence for a given (Z,m) mode is determined by the box operator D; im acting on it. This means 
that the asymptotic behavior of the R^ m mode is determined by the box operator, since the couplings are lost in that 
limit. Indeed, we find that solving the equation 

□i,m Ri, m = 0, (23) 

is enough to identify the main properties of the R^ m mode in regions far from the black hole. The numerical evolutions 
described in section |V] will show that solving this equation already gives a good description of the solution. This will 
be relevant in the study of the properties of late time behavior of the gravitational perturbation as tails (see @, EH 
for a review on the subject), or in the fast and efficient generation of waveforms. Notice that in the Schwarzschild 
case (i.e., a = 0) the modes decouple and we arc left with only the equation ((23)) for each mode. 

The system (|19p has another peculiarities. On one hand, since the gravitational waves only have quadrupole 
and higher contributions, the radial-temporal modes {i?i, m , Ro,m} are zero by construction and the first non-trivial 
equation would correspond to Z = 2, namely 

- a 2 A^ m d tt Ri.m - 4 a B 3 , m {^d tt - i <9 t ) R s>m + U 2 , m R 2 , m = . (24) 

On the other hand, the coupled system contains an infinity set of equations involving all possible values of I > 2. 
In order to solve it, the expansion (|14j) must be truncated at some order l m ax- From there on it is assumed that 
Rl>l ma!C m = 0, which prevents a overdeterminated system. The last equation of the series would look like: 



Htm a x,mRi ma x,m — 4:aBi max , m 4a- — —d tt - id t ) Ri max -i, m ~ a Ai maxtm d tt Rt maai -2,m = .(25) 

The modes Z = 3 and I = l max -i are also specials, since there still appear zero contributions. The rest of the modes 
in the system of evolution equations (|19|) are coupled by the next two values of the parameter I , above and below a 
given one. 



IV. THE FIRST ORDER EVOLUTION EQUATIONS 



In this section we will present two different formulations of the Teukolsky equations suitable for numerical evolutions. 
The first one is the first order reduction of the 1 + I coupled system (fl"9|) described in the previous section, where 
several Z— modes will be taken into account. In order to check the results, we will compare with a standard 2 + 1 
formulation to solve the gravitational perturbation equation. This second approach uses a different base of function 
to describe the angular dependence, leading to a two dimensional uncoupled system in the coordinates (r, 9). 
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A. The coupled one-dimensional system 

Since the Weyl scalar is complex, we first split the radial function into a real and imaginary part, that is, 

Rl,m = + Rl,m + i Rl,m ■ 

so that the eq. (|19|) decouples into the following second order system of equations: 



a M+i,m d t t Ri+2,m — 4a Bi + i^ n ( 4 a n^n du Ri+i, m ± d t T i?z+i,,„ ) + r D/, m Ri 



(26) 



(27) 



+2 a m (jffajjdt + d r -l) T Ri, m - 4 a B Lm (a c 
where we have defined the operator r dj m as 



<9t* ± Rl-l.m ± 3 t T -Ri_i. m ) — a 2 Al m <9 (t ± i?/_2.m — 0, 



r 2 + 2 Mr + — 
3 



1 + 2 



(l + 4) (i - 3) (I (I + 1) - 3m 2 ) \\ 9 



2 (2 r + M) — + 2 ( 2 r - M 



(21 + 3) {1 + 1)1(21-1) 
d n Mr + a 2 



, A V IMr 

9i 2 <9r 2 9i<9r 



(28) 



9r 



(Z-l) (Z + 2), 



The reduction from second order in space to first order can be achieved by defining the functions 



^l,m — dr ± Rl,m, 



d t ± Run+(3 ± ^l,m, 



(29) 
(30) 



For the Schwarzschild case ft is directly the shift vector and the above definition is useful to align the wave speed 
to the light cones. For the Kerr solution this is not the case but we explore the effect of varying it, the resultant 
wave form is insensitive to our choice and we set the easiest value of ft — 1, but we leave it as a free parameter in the 
equations. 

From these definitions we then obtain 



dt ± Ri, m — ± II/ i „ l — ft ^l } m, 

a t ± */, m - d r ( ± n i , m -/? ± * ; , m ), 

d u ± Ri, m = d^Hun-ftdr ( ± n ; , m -/? ± * ; , m ) 
Using this last equation we obtain the following first order equations: 



[(T 1 ,T 2 ,T 3 ,T i ,T 5 )d t + (c 1 ,c 3 ,c 5 ,c 7 ,c 9 ,)d r ] 



I ± n /+2 , m \ 

± n z , m 
± n / _ 1 . m 
V ± n l - 2 ,«, J 



[(c 2 ,C 4 ,C 6 ,C 8 ,Cio)9r + (Si, 5 3 , 57,512,514) 



(31) 
(32) 
(33) 



/ ±^1+2,™ \ 
± *i+l,n» ' 



+ (0,5 2 ,5 6 ,5n,0) 



T n z _i. m 

V ^l-2.m J 



+ (0,5 4 ,5 8 ,5i 3 ,0) 



/ T ^l+2,m \ 

T *Z-l.m 
V T *l-2,n J 



S 5 ±n ; , m + 5 9 ± R Um + 5i TR Um = 0, (34) 



which are the final first order equations. The coefficients are defined as 



T\ — — a 2 Ai + 2,m, T 2 — — 16 a 2 Bi + i, m jppn \ 1 ^3 

Ta = -16 a 2 B 



r 2 + 2 M r 



, , „ (Z+4) (1-3) (l(l+l)-3rn 2 ) 
1 + Z (2i+3) (l+l) I (21-1) 



l,m (; + i) (l-l) 1 T 5 — a 2 A^ m 

c l = -T 1 p, c 2 =T 1 ft 2 , c 3 = -T 2 ft, c± = T 2 ft 2 , c 5 =T 3 (4Mr-ft), c 6 = T 3 (ft 2 - 4 M r ft + A) 

c 7 = -T 4 ft, c s =T 4 ft 2 , c 9 = -T 5 /3, c 10 = T 5 ft 2 , 
S 1 =T 1 ftd r ft, S 2 = +AaB l+hm , S 3 =T 2 ftd r ft, Si = +AaB l+hm ft, 5 5 = 2 (2r + M) T 3 , 

S 6 = +8a 7J ^ JJ T 3 , 5 7 = ({ft-AMr) d r ft-2 (2r + M) ft + 2 (2r-M-£)) T 3 , 

_ ( nMr+a 



5 8 =T2am (^1-^^J T 3 , S 9 = (2^-((-l) (I + 2)jT 3 , 5 10 = ±2™fT : 
5 n =T4aB ; , m , S 12 =T A ftd r ft, 5 13 = ±4 a B, m ft, S u = T 5 ftd r ft. 
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By reducing the system to first order, we have enlarged our set of variables for each mode (I, to) from 2 (real and 
imaginary parts ± -Rz i7Ti ) to 6 ({^R^m, ± IIzm ^^im}). In order to solve this system of coupled equations we write it 
as a matrix equation, Mdt u-\-Ad r u-\-Bu — 0. The matrix M has to be inverted in order to obtain the evolution 
equation for each variable. 



B. The two-dimensional system 

We will follow here the approach given by Krivan et al.[3l but for the Tcukolsky function $ = r^^ instead of 
(r — ia cos #) 4 \E , 4^ 1 which will introduce some differences on the final equations. Assuming that the solution to eq. 
([9]) can be expressed in terms of the azimuthal modes, that is, 

$ = J2e im vS m (r,6,t), (35) 

rn 

it is obtained, for each azimuthal mode m, the following perturbation equation: 



d 2 


+ A 


d 2 


\-AMr 


d 2 


d 2 


dt 2 


dr 2 


dtdr 


de 2 + 


d 


cot t 




(M 


r + a 2 


. a to 

r 


dr 


de 







d_ 

dt 



(36) 



m 2 - 4 to cos 6 + 2 (l + cos 2 9) \ 

2 [ 2r-M +iam) — +cot9— +2 = i '- S m = 0. 

r J or oO \ r 2 r J sin 6 I 

The function S m is complex, so we separate it in its real and imaginary parts 

S m = S + + iS-, (37) 

and obtain two sets of equations 

P+S+ - P-S- = 0, 
P+S-+P-S+ = 0. 

where P + and P_ refers to the real and imaginary parts of the operator in equation (|36[) .As we did before, a first 
order reduction can be obtained by defining the variables Tl± as 

n± = d t S± + bd r S±, (38) 

where 

b = (-2Mr + J(2Mr) 2 + A(E + 2Mr)\ . (39) 

S + 2M r \ J 

In terms of the new variables n_|_,II_ the perturbation equation can be rewritten as: 

<9 t n± = ci5 r n± + c 2 d r S + ± c 3 d r S^ - c 4 n± ± c 5 U T - c 6 5± ± c 7 S T + O e ,<t>S±, (40) 
d t S± = U ± -bd r S ±7 (41) 

where now we have defined the following functions 



ci 



±— (2Mr + V((2Mr) 2 + A(E + 2Mr)) 
- 2M r V 



2(2r 2 - Mr - a 2 ) d 

c 2 = 6c 4 H — + c ^l~ b 

(lj + 2Mr)r dr 

2am 

C3 = -^2Mr- C " h 



Ci 
C6 



2(2r + M) 
~ S + 2Mr 

4a cos 
~ S + 2Mr 

r 2 m 2 csc 2 6< - 2(Mr - a 2 ) 2(1 + cos 2 6») - 4mcos( 



! (S + 2Mr) sin 2 0(E + 2Mr) 



« 1 n COte , 

" S + 2Mr^ + S + 2Mr < 
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Since we are using the Kerr-Schild coordinates, there is no need to introduce the tortoise coordinate to avoid sin- 
gularities, as our chart is regular at the apparent horizon. Our equation (|36p is similar in form at the one given by 
Campanelli et al. ||. 



V. NUMERICAL SIMULATIONS 



In this section we will evolve and analyze the gravitational response of the background spacetime for different initial 
perturbations. First, we will compute the truncated solution of the one-dimensional coupled system (|34[) for different 
Imax = 2,3 and 4, leading to a system with 6, 12 and 18 unknowns respectively These solutions will be compared 
with the solution from the usual two-dimensional perturbation equation (|4ip . This procedure will not only prove 
the validity of the Id description, but also will allow us to better understand the role of each mode as a function 
of the mode of the initial perturbation. Next, we will solve only for l max = 2 (i.e., the first order reduction of the 
harmonic equation (|23|1 ). comparing again with the solution from the two-dimensional system. As we will show, 
solving this single decoupled equation for each mode is an efficient way to describe quite accurately the gravitational 
perturbations. 

The numerical simulations have been performed with an extension of the code used in [B| for solving the perturbation 
equation in a Schwarzschild background. The equations are evolved by means of the Method Of Lines, with a third 
order Runge-Kutta for the time evolution and a fourth order differences to approximate the radial derivatives. A 
small amount of Kreiss-Oliger dissipation is added for numerical stability reasons. 

The reliability of the numerical implementation and the accuracy of the solutions was checked by performing several 
evolutions with different resolutions in order to determine the convergence of our results. We conclude that the solution 
shows a third order convergence, although higher resolutions are needed (i.e.,Ar ~ 0.04 for the finest grid) to keep 
this convergence order in the cases with high spins of the black hole (i.e., a = 0.999). 



A. Initial data with a 1=2 mode 



We study the response of a black hole with mass M = 1 to some initial infalling gravitational perturbation, for 
two different values of the rotational parameter a = 0.2 and a = 0.9. The angular dependence of the initial data 
is described by a single spin- weighted spherical harmonic mode, either with (I = 2, m = 0) or (I = 2, m = 2). The 
radial part has a Gaussian profile of the form g(r) = v4cxp(— (r — ro)/cr 2 ), with an amplitude, center and width 
respectively of A = 5 x 10 -3 , ro = 20M, a = 0.5. Its time derivative is assumed to be zero initially. This initial 
perturbation is evolved with the Id harmonic decomposition, solving the equation (|24|l with i?2,ro(f>^ = 0) = g(r). 
The one dimensional coupled system is evolved with l max = 2, 3 and 4 in order to calculate how fast the truncated 
solution converges to the complete one. 

In order to compare with the 2d approach, we have to write the corresponding initial data on this basis. For 
instance, in the case of the mode (I = 2, m = 0), that is 

S + (r,6,t = 0)=g(r)Y_ 2 20 (d), (42) 

with 

Y-2 2 \e) = l,/^sin 2 (0) . (43) 



We solve in a 2— dimensional grid the equation (|41j) . with the initial data given above (|42[) . In order to compare the 
solutions, we have to project the solution from the 2d system onto a spin weighted spherical harmonic basis, as long 
as we will compare the radial part of the wave form. In this way, we use the projection: 



2d R 2 , Q = j S + (r,d)Y w 2 (6) cM, (44) 

and then compare both signals, namely, the one produced by the Id harmonic decomposition and the one obtained 
by means of the standard 2d evolution. 

Our results for the (I — 2, m = 0) mode are shown in Fig.[TJ where it is displayed the signal measured by an observer 
located at r = 70M for different l max and the corresponding projection of the 2c? approach. There is a remarkable 
agreement in all the waveforms, showing that the results with l max — 2 are very similar to the ones with l max = 3 
and l max = 4 and to the 2d system. This means that the truncated series converge very fast and that the influence 
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FIG. 1: We present two values for the rotational parameter of the black hole a = 0.2 and a — 0.9 for an initial perturbation 
described by the I — 2, n er located at r — 70M, 

in the second the same s ained for a — 0.9 




100 200 300 

time 



FIG. 2: As in the previous figure, but for an initial perturbation described by the I = 2, m = 2 mode 

of the coupled neighbor modes is pretty small. Since solving equation ([24]) is much less expensive than solving the 2d 
code, this may allow us to use a finer grid with more resolution to describe the waveform. 

Our next example is to consider a initial perturbation with the same radial profile but with an angular dependence 
given by the (I = 2,m = 2) spin-weighted spherical harmonic mode. Notice that this is an extreme case since the 
coupling of modes depends strongly on m. The results are displayed in figure ([2]). In this case, there is a noticeable 
difference among the results obtained in the Id harmonic decomposition, depending on the number of variables we 
evolve, that is, on the modes that we let to be awaken. The contribution of the higher modes is evident and, as we 
increase number of modes that we let to be awaken, we get closer to the signal of the 2d approach. 

The results presented in this subsection not only indicate that the truncated solution of the Id coupled system 
describes accurately the waveforms of a gravitational perturbation, but also allow us to conjecture that the coupling 
among the modes is mainly determined not so much by the value of the rotation parameter of the black hole, but by 
the value of the azimuthal mode number m. 
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FIG. 3: It is shown in first panel the signal produced by the axial code, measured by an observer located at r — 80M in 
the equatorial plane 9 = 7r/2, and the rotational parameter of the black hole is a = 0.999. In the panel below it is shown a 
superposition of the projection of that signal in the Y20 and Y30 mode and the ones obtained by the Id harmonic decomposition. 

B. Initial data with two modes with different 1 

Next we will consider an initial perturbation whose angular dependence is described by a combination of the two 
modes I = 2 and I — 3. The problem is made more challenging by increasing the black hole rotational parameter to 
a = 0.999. In this case we evolve only within the Id harmonic decomposition (|23p for each mode separately. Since 
any perturbation can be decomposed in spin weighted spherical harmonics, this example gives us some insight about 
the general behavior on the interaction of modes initially awaken. The initial data is given by the two modes describe 
by the same radial function R2.0 = R3.0 = Acxp(— (r — ro)/a 2 ), with the time derivative of both modes set to zero 
initially. 

For the 2d approach we have the linear superposition of those modes 

S+(x,9,t = 0) = g(r) (F_ 2 2 » + F- 2 3 ») . (45) 

As before, we make the projection of the final signal in order to compare both approaches. The results arc plotted in 
Fig. ©. Notice that, despite the fact that two modes are initially present, we still are able to describe correctly their 
evolution by means of the Id harmonic decomposition, showing that it could be used instead of the usual procedure 
of solving the two dimensional description for S m , even for very large value of the rotational parameter of the black 
hole. 



VI. CONCLUSIONS 

Wc have described the perturbation equation for a Kerr black hole in Kcrr-Schild coordinates. By finding a suitable 
tetrad, we have been able to obtain a perturbation equation the Weyl scalar $4, which is amiable to be described 
by the spin weighted spherical harmonics instead of the spheroidal ones used in other previous works. This fact 
allowed us to use the well known properties of the spin weighted spherical harmonics for separating the angular 
dependence from the temporal-radial one. Finally, we obtained a system of equations on the coupled coefficients of 
the perturbation equation. This coupling is given by the closest modes, two above and two below (in terms of I), and 
there is not an explicit mixture with respect to the mode number m. Achieved in this way, a separable description 
for the gravitational perturbation in the Kerr background, in the time domain. 

Once we obtained this system of equations, we proceed to test its actual behavior by making several numerical 
codes. We defined new variables that allowed us to re-write the coupled evolution equation as a first order system. 
We solve it numerically by splitting the real and imaginary parts. We compare the results obtained by this method 
with the standard one, writing the perturbation equation by a decomposition in term of the axial modes only. When 
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we compare our results with the coupled modes with the axial one, we found that the coupling manifest itself in an 
stronger way by the mode number m, rather than by the values of the angular parameter of the black hole. Indeed, 
in the examples that we considered, we saw that the wave form obtained by solving simply the evolution equation for 
the corresponding mode without considering the influences of the neighborhood modes coincides with great accuracy 
with the wave form obtained with the axial evolution. We also showed that the QNM frequencies of the ring down 
phase of the wave form, agree very well with the known results. 

In this way we have found a new, suitable approach that allows to generate accurate wave forms in a Kerr background 
in efficient way. 

Acknowledgments 

DN acknowledges DAAD and DGAPA-UNAM grants for partial support, as well to Lidia Rojas for her support 
during the elaboration of the present work. JCD acknowledges CONACyT support. 

Appendix A: Source terms 

The explicit form of the operators acting on the perturbed stress energy tensor. 

T = T kk T^ kk + T km ' T^\ m > + T m * m * T« m , m , , (Al) 

with 

T kk = -(<T +7712 (3a + /T +4tt-t*)) (<5*+77 12 (2a + 2/T-r*)), 
T fcm * = (A + m2 (4 M + M *+3 7 - 7 *)) (S* + 2 m2 (a - r*)) + (5* + jji 2 (3a + /3* + Air - r*)) {A + 2 Vl2 (^+7)), 
T m ' m * = — (A + m (4 n + li* + 3 7 - 7*)) ( A + ma (n* + 2 7 - 2 7 *)) . (A2) 

We recall that the T^\b are the projections of the perturbed source term, T' 1 ^,, on the null tetrad. 

In the Kerr background, with our choice of coordinates and null tetrad, the operators of T, defined by Eqs. (|A1IA2[) . 
have the following explicit form: 

s 2 



rrK k _ _ 

~ 2 
r a 2 



T 



d 2 



i) 



B.^o-a' sin'e \ —-2{2^ + ^) \-~n*\\-2iasin9 I — — 2 /i — /i' , 



d 



ft. - — 

dt dr 

d 2 



d_ 

dt 



d_ 

dr 



— 2 

dt 2 dtdr dr 2 



d 2 ( d d 



(A3) 



which trivially reduce to the corresponding expressions in the Schwarzschild case, [f|. 



Appendix B: Scalar field 

As mentioned in the text, a similar procedure can be use to deal with the motion equation for a massless scalar 
field, which satisfies the corresponding Klein-Gordon equation, V Q V a $ = AirT. Indeed, in the Kerr background, 
given by the line element Eq. Q, the Klein-Gordon equation takes the explicit form: 



P + DflJ = 4ttT. 



where we have defined <p = r$, and the operators: 



□ = - (r 2 +a 2 cos 2 1 



d 2 d 2 d 2 
' dt 2 dr 2 dtdr 



~ 2M Ft +2 ( M -V)d-r + 



a,, 



d 2 



1 



d 2 



cot I 



d_ 

dr 
d 



M r - a 2 



a d 



d9 2 sm 2 9 dip 2 ' """" d9 r dip' 



(Bl) 



2a 



drdip 



(B2) 
(B3) 
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The form of the angular operator is very suggestive to expand the function in terms of the spherical harmonics, 

= X^i,m(t,r)ro ,,m (M> ( B4 ) 



In 



and as, Dg v Y l ' m = (-1 (I + 1) + 2 i^f) Y l ' m , and 2 a ^ Y Q l ' m = Y l ' m 2 i a . We get again an almost radial 
temporal equation with only an angular coefficient a 2 cos 2 9 of the second temporal derivative. 

Using the orthonormalization properties of the spherical harmonics and the fact that § Yo Lm Yq 1 ' m Yq 20 = 





V 




(1 V 2 \ 


i 





! 


I m to' / 



where we are using the Wigner-3 matrices |15l | , again we obtain 
that we can completely get rid of the angular dependence in favor of a system with coupled modes: 

d 2 d 2 

n<j) = —a Asi +2 ,m-Qp El+2,m + □« El, m - a Asi. m — El-2,m (B5) 

where we have now defined 



(I + m)(l + to - l)(l - m)(l - to - 1) 
M, m = — K \ (2Z + l)(2Z-3) (B6) 




l 2 + l-'im 2 \\ d 2 A d 2 d 2 
^ + A—r +4:Mr- 



(21 + 3) (21 -I) J J dt 2 dr 2 dtdr 

„ . r 9 „ /, , a 2 \ 3 n M r — a 2 — imar ,,, 

-2M— + 2 M + ia?7i — + 2 = i (I + 1) . 

ot \ r J or r z 

In this way, we proved that the Klein-Gordon equation for a massless scalar field in a Kerr background can be 
described, as the gravitational perturbation, by a radial-temporal system of equations for the coupled modes. This 
system of equations has similar properties to the ones mentioned with respect to the gravitational perturbation, and 
is remarkable that for the scalar case the coupling of the modes involve only two neighbor modes, the second neighbor 
above and below a given 1 mode number. 

The description of the evolution of the scalar field presented here, clarifies some aspects regarding the coupling 
among modes, which have been discussed in the literature, see for instance [H, and will be useful to study the 
properties of the scalar field evolving in this background, such as the late time behavior known as tails. 



Appendix C: Electromagnetic field 

The description of the electromagnetic field in a Kerr background, can also be described within the present for- 
mulation. Starting from the Maxwell's equations, F^ v . v = AttJ^, and considering that the Faraday tensor can be 
expressed in terms of the tetrad vectors as [l6| : 



(fam*^k u] + fal lfl m^ +<t>i (k^l v] + TO^TO* iy] )) +c.c, (CI) 



where fa = F^ v to", fa = F liv ra*^k v and fa — \F }1V k v + m*^ m u ) arc complex scalars. Projecting the 
Maxwell's equations onto a null tetrad, we obtain for type D spaces, the following set of equations for the complex 
scalars: 

faD-2p) <f>i-(rj8* -2a + Tr) fa = 4tt J u (C2) 

(r)5 + 2p-T) fa-(r)A + 2p,) fa = 4tt J kl (C3) 

(r]S-2T) fa-(r]A + fx-2y) fa = 4ttJ to , (C4) 

(77D -p + 2e) fa - (r)6* + 2tt) fa = 4nJ m «, (C5) 

where J a are the projections of the current vector, J M on the respective null vector, and 77 = ±1 depending on the 
choice of signature (for the present work, 77 = —1). 

Following the usual procedure, [3], of acting with specific operators on the equations and using the commutation 
relations, we get the following equations for fa and fa '■ 

[( V r>-2p-p*-e + e*) (nA + (i-2i)-( V 6-a* -P + tt* -2t) (nS* -2a + n)] fa = 4ttJ , (C6) 
[(77 A + 2 /j + M * +7-7*) faD -p + 2e) - (77 J* + a + (3* + 2tt-t*) (r?5 + 2/3-r)] fa = 4ttJ 2 , (C7) 
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and wc have defined 

Jo = (r]5-a*-p + n*-2T) J, - (r?D - 2p - p* - e + e*) J m , (C8) 
J 2 = ( J? A + 2 /1 + M , +7-7 , )4-(^* + " + /5 , +2^-t*)4 (C9) 

A direct substitution of the operators and spinor coefficients for the Kerr background described by the line element 
given in Eq. ([3]) using the null tetrad given by Eq. (j4]), after simplification produces the following equations: 



UiD + iiD^] <j>±i =4ttJ±i. (CIO) 

have defined 4>\ = (r — i a cos(6* 

operators: 



where we have defined $± = (r — i a cos(9)) 2 4>q, (/>_i — r 4> 2 , J\ = —2 £ (r — ia cos(6*)) 2 J , J_i = —2 r £ J 2l and the 



i Q = -(r 2 + a 2 cos 2 9 + 2Mr) |L + A |1 + 4Mr ^ + 2 a ^ - 2 (r + i a co S 0) |, (Cll) 

^2 ^2 ^2 ^2 ^ 

-iD = - (r 2 + a 2 cos 2 9 + 2 Mr) —~ + A + 4 M r — — + 2 a — — + 2 [r + ia cos (9) — + 
v ' ot ar 2 utor drop) at 

( r ) <9r r dip r 2 ' ^ ^ 



<9 2 i a 2 „ 9 , n .cos9 a i 

<96» 2 sin 2 dip 2 39 sm 2 sm 2 ° 



Thanks to the tetrad we are using, Eq. (U), we again obtain that the main angular operators, ±ide v , are such that 
the spin weighted spherical harmonics, Y±i ,m , are their eigenfunctions. Indeed, ±i\Dg v Y±i" m = —I (I + 1) Y±i' m . 
As K)-!''" 1 = im Y±i' m , we expand the electromagnetic scalars in terms of such harmonics: 

hi = E ±i^,«(*. r ) ^±i' ,m (^ p)- (C14) 

Substituting in the field equations, Eqs. ([C10[) we get, as for the other fields, equations where the angular dependence is 
reduced to the coefficient cos 2 6 for ±1 7?.; m , and the coefficient cos0 for ±iR-i, m - We expand these trigonometric 
functions as spherical harmonics and, using the orthonormal properties with the Wigner matrices, we are able to 
express this final angular dependence in favor of a coupling of the modes, obtaining the following system of equations: 

,m + ^lO^m zpllZl^n (C15) 



where wc have now defined 



2 _ _L / (' + ^Kj + m- lKj - m)(Z - m - 1)(Z + 1)(Z - 2) 
*' m " (2Z-1)V (2Z + 1) Z (Z - 1) (2Z - 3) ' 1 J 



5 TO / (Z + m)(l-m) 

- T V (2Z + 1)(2Z-1)' (C17) 

- / 2 a 2 A , (« (/ + 1) -3)(f(f + l) -3m 2 )\\ d 2 A 5 2 , „, 8 2 

&l, m = - r 2 + 2Mr + - 1 + 4^ > ■ >\ ±- > ^ + A ^+ 4M ^ 



3 V (21 + 3) (21 - 1) J J dt 2 dr 2 dtdr 

i 2 „„ a2 A < (Z(Z + l)-3)(Z(Z + l)-37n 2 )\\ d 2 A 3 2 4 „ 9 2 
= - |r +2Mr + i ( 1 + 4 (2Z + 3)(2i-l) "J J ^ + A ^ +4M ^ + 

2 (r + 4ia— ^— ^) |- + 2 (r- — +iam) - 2 im- + 2% - I (I + 1) . (C18) 
1(1 + 1)/ oi \ r J or r r 2 
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Thus, the dynamics of the electromagnetic field in a Kerr background can be described, as the other fields, by a 
radial-temporal system of equations for the coupled modes. This system of equations has similar properties to the 
ones mentioned with respect to the gravitational perturbation. 



[1] S. A. Teukolsky. Perturbations of a rotating black hole i. fundamental equations for gravitational, electromagnetic, and 

neutrino-field perturbations. Astroph. J., 185:635-647, 1973. 
[2] R. Price. Nonspherical perturbations of relativistic gravitational collapse. I. Scalar and gravitational perturtbations. Phys. 

Rev. D, 5:2419, 1972. 

[3] R. J. Gleiser, R. H. Price, and J. Pullin. Late-time tails in the kerr spacetime. Class. Quantum Grav., 25:072001, 2008. 
[4] L. M. Burko and G. Khanna. Late-time kerr tails revisited. Class. Quantum Grav., 26:015014, 2009. 

[5] J. C.Degollado, D. Nunez, and C. Palenzuela. Signatures of the sources in the gravitational waves of a perturbed 

schwarzschild black hole. General Relativity and Gravitation, To appear, 2009. 
[6] C. Moreno and D. Nunez. Gravitational perturbations of the kerr black hole due to arbitrary sources. Int. J. of Mod. 

Physics, D 11:1331-1346, 2001. 

[7] A. Zenginoglu, D. Nunez, and Sasha Husa. Gravitational perturbations of schwarzschild spacetime at null infinity and the 

hyperboloidal initial value problem. Class. Quantum Grav., 26:035009, 2009. 
[8] M. Campanelli, G. Khanna, P. Laguna, J. Pullin, and M. P. Ryan. Perturbations of the kerr spacetime in horizon- 
penetrating coordinates. Class. Quantum Grav., 18:1543-1554, 2001. 
[9] E. Newman and R. Penrose. Note on the bondi-metzner-sachs group. J. Math Phys., 7:863-870, 1966. 
[10] J. N. Goldberg, A. J. Macfarlane, E. Newman, F. Rohrlich, and E. C. G. Sudarshan. Spin s spherical harmonics and eth. 

J. Math Phys., 8:2155-2161, 1967. 
[11] E. W. Leaver. An analytic representation for the quasi-normal modes of kerr balck holes. Proc. R. Soc. Lond., A 402:285- 
298, 1985. 

[12] M. Ruiz, M. Alcubierre, D. Nunez, and R. Takahashi. Multipole expansion for energy and momenta carried by gravitational 

waves. Gen. Relativ. and Gravit., 40:1705-1729, 2008. 
[13] A. Zenginoglu. A hyperboloidal study of tail decay rates for scalar and yang-mills fields. arXiv:0803.2018vl, gr-qc, 2008. 
[14] W. Krivan, P. Laguna, Ph. Papadopoulos, and N. Andersson. Dynamics of perturbations of rotating black holes. Phys. 

Rev., D 56:3395-3404, 1997. 
[15] M. Alcubierre. Introduction to 3 + 1 numerical relativity. Claredon Press, Oxford, 2007. 
[16] S. Chandrasekahar. The Mathematical theory of black holes. Clarendon Press Oxford, Oxford, UK, 1983. 



